%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%THIS FILE CONSTRUCTS SOME ESTIMATED PARAMETERS FROM THE MODEL IN 
%"Evidence of Neighborhood Effects from Moving to Opportunity: 
%     LATEs of Neighborhood Quality"
%by Dionissi Aliprantis and Francisca G.-C. Richter
%SEE THE INCLUDED readme.txt FILE FOR AN EXPLANATION OF THIS FILE 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%IMPORTING ESTIMATION RESULTS AND DATA
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
load('xhat.mat');
x=xhat;
    data=importdata('C:\MTO_1\matlab\matlab_data.txt');

N=length(data(:,15));              %NUMBER OF OBSERVATIONS
T=2;                        %NUMBER OF TYPES (t=1 for S8, t=2 for Experimental)

cutpoint=zeros(N,11);
cutpoint(:,2)=data(1,15);
cutpoint(:,3)=data(2,15);
cutpoint(:,4)=data(3,15);
cutpoint(:,5)=data(4,15);
cutpoint(:,6)=data(5,15);
cutpoint(:,7)=data(6,15);
cutpoint(:,8)=data(7,15);
cutpoint(:,9)=data(8,15);
cutpoint(:,10)=data(9,15);
cutpoint(:,11)=50;

data=data(:,1:14);
data( any(data==-9,2), : ) = [];



%%%%%%%%%%%%%%%%
%LABEL VARIABLES
%%%%%%%%%%%%%%%%
D10         =data(:,1);     %D \in {1, 2, ..., 10}
HHnoteens   =data(:,2);     %X_1 
HHnofam     =data(:,3);     %X_2 
HHvictim    =data(:,4);     %X_3 
basequal    =data(:,5);     %X_4 
sitedummy2  =data(:,6);     %BOSTON  SITE DUMMY
sitedummy3  =data(:,7);     %CHICAGO SITE DUMMY
sitedummy4  =data(:,8);     %LA      SITE DUMMY
sitedummy5  =data(:,9);     %NYC     SITE DUMMY
exp         =data(:,10);    %INDICATOR FOR MTO EXPERIMENTAL VOUCHER 
S8          =data(:,11);    %INDICATOR FOR SECTION 8 VOUCHER
weight      =data(:,12);    %WEIGHTS
mover       =data(:,13);    %INDICATOR FOR MOVE WITH VOUCHER
qstar       =data(:,14);

N=length(D10);              %NUMBER OF OBSERVATIONS

cutpoint=cutpoint(1:N,:);




%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%CREATING LATENT INDEX COMPONENT BY COMPONENT
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%
%mu(X)
%%%%%%
bHHnoteens   =x(1);
bHHnofam     =x(2);
bHHvictim    =x(3);
bbasequal    =x(4);
bsitedummy2  =x(5);
bsitedummy3  =x(6);
bsitedummy4  =x(7);
bsitedummy5  =x(8);

mu=zeros(N,1);
mu(:)=bHHnoteens*HHnoteens(:) +bHHnofam*HHnofam(:) +bHHvictim*HHvictim(:) +bbasequal*basequal(:) ...
    +bsitedummy2*sitedummy2(:) +bsitedummy3*sitedummy3(:) +bsitedummy4*sitedummy4(:) +bsitedummy5*sitedummy5;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%COST FUNCTION (ie, CUTOFFS)
%
% !!!!!! OJO !!!!!!!!!!!!!
%
%HERE IT IS COST OF GOING FROM j-1 TO j
%IN PAPER IT IS COST OF GOING FROM j TO j+1
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
bcutpoint(1)   =-Inf;
bcutpoint(2)   =x(9);
bcutpoint(3)   =x(10);
bcutpoint(4)   =x(11);
bcutpoint(5)   =x(12);
bcutpoint(6)   =x(13);
bcutpoint(7)   =x(14);
bcutpoint(8)   =x(15);
bcutpoint(9)   =x(16);
bcutpoint(10)  =x(17);

C=zeros(N,11);               %COST CUTOFFS
%should be index=-Inf --> C(:,1)=-Inf;
for j=1:10
    C(:,j)=bcutpoint(j);
end
C(:,11)=Inf;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%COST-REDUCTION FROM SECTION 8 VOUCHERS (FOR HIGH-TYPES)
%
% !!!!!! OJO !!!!!!!!!!!!!
%
%HERE IT IS COST-REDUCTION OF GOING FROM j-1 TO j
%IN PAPER IT IS COST-REDUCTION  OF GOING FROM j TO j+1
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%g8cutpoint(1)    =x(18); - IRRELEVENT WITH NEW j INDEXING
g8cutpoint(2)    =x(18);
g8cutpoint(3)    =x(19);
g8cutpoint(4)    =x(20);
g8cutpoint(5)    =x(21);
g8cutpoint(6)    =x(22);
g8cutpoint(7)    =x(23);
g8cutpoint(8)    =x(24);
g8cutpoint(9)    =x(25);
g8cutpoint(10)   =x(26);

g8sitedummy2     =x(27);
g8sitedummy3     =x(28);
g8sitedummy4     =x(29);
g8sitedummy5     =x(30);

gamma8=zeros(N,11);          %EFFECT OF VOUCHER INCLUDING UNOBSERVED HETEROGENEITY
for j=2:10
    gamma8(:,j)=g8cutpoint(j) +g8sitedummy2.*sitedummy2(:) +g8sitedummy3.*sitedummy3(:) ...
                                           +g8sitedummy4.*sitedummy4(:) +g8sitedummy5.*sitedummy5(:);
end
gamma8(:,11)=0;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%COST-REDUCTION FROM EXPERIMENTAL VOUCHERS (FOR HIGH-TYPES)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%g1cutpoint(1)   =x(31); - IRRELEVENT WITH NEW j INDEXING
g1cutpoint(2)   =x(31);
g1cutpoint(3)   =x(32);
g1cutpoint(4)   =x(33);
g1cutpoint(5)   =x(34);
g1cutpoint(6)   =x(35);
g1cutpoint(7)   =x(36);
g1cutpoint(8)   =x(37);
g1cutpoint(9)   =x(38);
g1cutpoint(10)  =x(39);

g1sitedummy2  =x(40);
g1sitedummy3  =x(41);
g1sitedummy4  =x(42);
g1sitedummy5  =x(43);

gamma1=zeros(N,11);          %EFFECT OF VOUCHER INCLUDING UNOBSERVED HETEROGENEITY
for j=2:10
    gamma1(:,j)=g1cutpoint(j) +g1sitedummy2.*sitedummy2(:) +g1sitedummy3.*sitedummy3(:) ...
                                           +g1sitedummy4.*sitedummy4(:) +g1sitedummy5.*sitedummy5(:);
end
gamma1(:,11)=0;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%MARGINAL BENEFIT (ie, LATENT INDEX) OF INDIVIDUAL n FROM MOVING FROM LEVEL j-1 TO LEVEL j 
%     (mu + gamma8_j^H*Z8 + gamma1_j^H*Z1 - C_j - V)
%
% !!!!!!!!! OJO !!!!!!!!!!!!!!!!!!!
%
%IN PAPER IT IS
%MARGINAL BENEFIT (ie, LATENT INDEX) OF INDIVIDUAL n FROM MOVING FROM LEVEL j TO LEVEL j+1 
%     (mu + gamma8_j^H*Z8 + gamma1_j^H*Z1 - C_j - V)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
index=zeros(N,10);          %LATENT INDEX FOR NBD QUALITY SELECTION
index(:,1)=-C(:,1);
for j=2:11
    index(:,j)=mu(:)+exp(:).*mover(:).*gamma1(:,j)+S8(:).*mover(:).*gamma8(:,j)-C(:,j);
end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%CREATING TYPE PROBABILITIES
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%
%SECTION 8
%%%%%%%%%%
b8HHnoteens   =x(44);
b8HHnofam     =x(45);
b8HHvictim    =x(46);
b8basequal    =x(47);
b8sitedummy2  =x(48);
b8sitedummy3  =x(49);
b8sitedummy4  =x(50);
b8sitedummy5  =x(51);

b8constant    =x(52);

index8=zeros(N,1);          %LATENT INDEX FOR SECTION 8 TYPE PROBABILITY
index8(:)=b8constant+b8HHnoteens*HHnoteens(:) +b8HHnofam*HHnofam(:) +b8HHvictim*HHvictim(:) +b8basequal*basequal(:) ...
           +b8sitedummy2*sitedummy2(:) +b8sitedummy3*sitedummy3(:) +b8sitedummy4*sitedummy4(:) +b8sitedummy5*sitedummy5(:);


%%%%%%%%%%%%%
%EXPERIMENTAL
%%%%%%%%%%%%%
b1HHnoteens   =x(53);
b1HHnofam     =x(54);
b1HHvictim    =x(55);
b1basequal    =x(56);
b1sitedummy2  =x(57);
b1sitedummy3  =x(58);
b1sitedummy4  =x(59);
b1sitedummy5  =x(60);

b1constant    =x(61);

index1=zeros(N,1);          %LATENT INDEX FOR EXPERIMENTAL TYPE PROBABILITY
index1(:)=b1constant+b1HHnoteens*HHnoteens(:) +b1HHnofam*HHnofam(:) +b1HHvictim*HHvictim(:) +b1basequal*basequal(:) ...
           +b1sitedummy2*sitedummy2(:) +b1sitedummy3*sitedummy3(:) +b1sitedummy4*sitedummy4(:) +b1sitedummy5*sitedummy5(:);

%%%%%%%%%%%%%%%%%%%%%%%%%%
%COVARIANCE MATRIX "SIGMA"
%%%%%%%%%%%%%%%%%%%%%%%%%%
sigmaDD=1;
sigma88=1;
sigma11=1;

sigmaD8=x(62);
sigmaD1=x(63);

SIGMA8=zeros(2,2);
SIGMA8=[sigmaDD sigmaD8 ; 
        sigmaD8 sigma88  ; ];

SIGMA1=zeros(2,2);
SIGMA1=[sigmaDD sigmaD1 ; 
        sigmaD1 sigma11  ; ];

    
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%POST-ESTIMATION CONSTRUCTIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    Cq=zeros(N,1);
    Cq(:)=Cq(:)+(-2+(qstar(:)-cutpoint(:,1)).*(C(:,1+1)+2)./(cutpoint(:,1+1)-cutpoint(:,1))).*(D10(:)==1);
for j=2:9
    Cq(:)=Cq(:)+(C(:,j)+(qstar(:)-cutpoint(:,j)).*(C(:,j+1)-C(:,j))./(cutpoint(:,j+1)-cutpoint(:,j))).*(D10(:)==j);
end
    Cq(:)=Cq(:)+(C(:,10)+(qstar(:)-cutpoint(:,10)).*(3-C(:,10))./(cutpoint(:,10+1)-cutpoint(:,10))).*(D10(:)==10);
scatter(qstar(:),Cq(:))

    gamma1q=zeros(N,1);
for j=2:9
    gamma1q(:)=gamma1q(:)+(gamma1(:,j)+(qstar(:)-cutpoint(:,j)).*(gamma1(:,j+1)-gamma1(:,j))./(cutpoint(:,j+1)-cutpoint(:,j))).*(D10(:)==j);
end
    gamma1q(:)=gamma1q(:)+(gamma1(:,10)+(qstar(:)-cutpoint(:,10)).*(gamma1(:,10)-.2-gamma1(:,10))./(cutpoint(:,10+1)-cutpoint(:,10))).*(D10(:)==10);
    %OJO - FOR NOW JUST KEEPING CONSTANT ACROSS D==1    
gamma1q(:)=gamma1q(:)+gamma1(:,2).*(D10(:)==1);
%    gamma1q(:)=gamma1q(:)+(gamma1(:,2)+(qstar(:)-cutpoint(:,2)).*(gamma1(:,3)-gamma1(:,2))./(cutpoint(:,3)-cutpoint(:,2))).*(D10(:)==1);

gamma8q=zeros(N,1);
%ORIGINAL
%     gamma8q(:)=gamma8q(:)+gamma8(:,2)+(qstar(:)-cutpoint(:,2)).*(gamma8(:,3)-gamma8(:,2))./(cutpoint(:,3)-cutpoint(:,2))).*(D10(:)==1);
% for j=2:9
%     gamma8q(:)=gamma8q(:)+(gamma8(:,j)+(qstar(:)-cutpoint(:,j)).*(gamma8(:,j+1)-gamma8(:,j))./(cutpoint(:,j+1)-cutpoint(:,j))).*(D10(:)==j);
% end
%     gamma8q(:)=gamma8q(:)+(gamma8(:,10)+(qstar(:)-cutpoint(:,10)).*(gamma8(:,10)-.2-gamma8(:,10))./(cutpoint(:,10+1)-cutpoint(:,10))).*(D10(:)==10);
for j=2:9
    gamma8q(:)=gamma8q(:)+(gamma8(:,j)+(qstar(:)-cutpoint(:,j)).*(gamma8(:,j+1)-gamma8(:,j))./(cutpoint(:,j+1)-cutpoint(:,j))).*(D10(:)==j);
end
    gamma8q(:)=gamma8q(:)+(gamma8(:,10)+(qstar(:)-cutpoint(:,10)).*(gamma8(:,10)-.2-gamma8(:,10))./(cutpoint(:,10+1)-cutpoint(:,10))).*(D10(:)==10);
%OJO - FLIPPED AT END - cutpoint(:,2) > qstar(:) when D==1
%KEEPING SLOPE gamma8(:,3)-gamma8(:,2) FROM D==2 
%    gamma8q(:)=gamma8q(:)+gamma8(:,2)+((cutpoint(:,2)-qstar(:)).*(gamma8(:,3)-gamma8(:,2))./cutpoint(:,2)).*(D10(:)==1);
%OJO - FOR NOW JUST KEEPING CONSTANT ACROSS D==1
gamma8q(:)=gamma8q(:)+gamma8(:,2).*(D10(:)==1);
    %+((cutpoint(:,2)-qstar(:)).*(gamma8(:,3)-gamma8(:,2))./cutpoint(:,2)).*(D10(:)==1); 

V=zeros(N,1);
V(:)=mu(:)+gamma8q(:).*(S8(:)==1).*mover(:)+gamma1q(:).*(exp(:)==1).*mover(:)-Cq(:);

 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%SAVING RESULTS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
site=zeros(N,1);
site(:)=1+1.*(sitedummy2(:))+2.*(sitedummy3(:))+3.*(sitedummy4(:))+4.*(sitedummy5(:));
    
results=zeros(N,20);
for i=1:62
results(i,1)=x(i);
end
results(:,2)    =mu(:);
results(:,3)    =qstar(:);
results(:,4)    =Cq(:);
results(:,5)    =gamma1q(:);
results(:,6)    =gamma8q(:);
results(:,7)    =site(:);
results(:,8)    =V(:);
results(:,9)    =weight(:);
results(:,10)   =basequal(:);
results(:,11)   =exp(:);
results(:,12)   =S8(:);
results(:,13)   =mover(:);

csvwrite('estimation_results.txt', results)






